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The ground state of the square lattice bilayer quantum antiferromagnet with nearest (Ji) and 
next-nearest (J2) neighbour intralayer interaction is studied by means of the dimer expansion 
method up to the 6-th order in the interlayer exchange coupling J3. The phase boundary between 
the spin-gap phase and the magnetically ordered phase is determined from the poles of the biased 
Pade approximants for the susceptibility and the inverse energy gap assuming the universality 
class of the 3-dimensional classical Heisenberg model. For weak frustration, the critical interlayer 
coupling decreases linearly with a(= J2/J1). The spin-gap phase persists down to J3 = (single 
layer limit) for CWS^oKO.GS. The crossover of the short range order within the disordered phase 
is also discussed. 

KEYWORDS: bilayer Heisenberg antiferromagnet, frustration, Pade approximant, dimer expansion method, spin- 
gap state 
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§1. Introduction 

The spin-1/2 square lattice Heisenberg model is 
now widely believed to have an antifcrxamagnetic long 
range order in the ground stateHrcra'cP It is, how- 
ever, expected that the strong quantum fluctuation in 
this system may lead to the destruction of the long 
range order with the help of some additional mech- 
anism. In this context, the square lattice antifer- 
romagnetic Heisenberg model with nearest and next- 
nearesL £xchang£i r b±Hr^tifjiL ^pmafter called Ji-J 2 
moddffBlWBBBBeBa bi- 

layer Heisenberg modeOBBBBBBBE] H have 
been studied extensively. Considering the difference of 
the nature of the mechanism leading to the spin-gap 
phase in these two models, it must be most interesting 
to study their interplay in the bilayer J\-J<i modelB3'c3 

In the bilayer model, if the interlayer antiferromag- 
netic coupling is strong enough, the spins on both layers 
form interlayer singlet pairs and the quantum fluctua- 
tion is enhanced leading to the quantum disordered state. 
The dimc.r_c:y n a n si op-study of this model has been quite 
successfuH'BBESEiP and it is shown that the transi- 
tion between the Neel phase and the spin-gap phase be- 
longs to the universality class of 3-dimcnsional classical 
Heisenberg model. This resulLis also confirmed by quan- 
tum Monte Carlo simulation.^ 

On the other hand, in the J1-J2 model, the com- 
petition between the nearest neighbour interaction J\ 
and the nearest neighbour interaction J 2 introduces the 
frustration in the spin configuration which enhances the 
quantum fluctuation. The conclusion about the pres- 
ence of the quantum disordered state in this model is, 
however, still controversial even in the most frustrated 



regime. 

In order to apply the dimer expansion method to the 
single layer J1-J2 model, it is inevitable to start with the 
dimer configurations which break the txanslational sym- 
metry as an unperturbed ground state.utr In the bilayer 
J1-J2 model, the unperturbed ground state can be taken 
as the interlayer dimers and the translational symmetry 
of the original Hamiltonian is preserved throughout the 
calculation. Therefore the bilayer model is more suit- 
able for the dimer expansion study than the single layer 
model. It is also possible to get insight into the phase 
transitions in the single layer model by investigating of 
the asymptotic behavior in the limit of vanishing inter- 
layer interaction. 

This paper is organized as follows: The bilayer J1-J2 
model Hamiltonian is introduced ii^tijCjaext section. In 
§3, the dimer expansion methodol§c2rti3 is applied to 
this model and the phase diagram is determined using 
the biased Pade analysis. The last section is devoted to 
summary and discussion. 

§2. Bilayer Ji-J 2 Model 

The Hamiltonian of the bilayer Ji-J 2 model is given 
as follows, 



H = J t £ (Sfsf + Sfsf) 

<i,j>nn 

+ j 2 J2 (sfsf + sfsf) 

<j,i>„»n 
+ J 3 J2sfsf, (2.1) 



where S i and Si are the spin operators with magnitude 
1/2 on the i-th site of the layer A and B, respectively. 
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The expression > and 



denote the summa- 



E 

<»ji>nn <l,j>n-nn 

tion over the intralayer nearest neighbour pairs and next 
nearest neighbour pairs, respectively. The last term rep- 
resents the interlayer coupling. All exchange couplings 
are assumed to be antiferromagnetic. In the following, 
we denote the ratios J2/J1 — cv. and J3/J1 = (3. In the 
classical limit, the ground state is the Neel state or the 



collincar 
a c = 0.5 



te according as a < a c or a > a c where 



§3. Dimer Expansion Method 

In the absence of the intralayer coupling, the ground 
state is the assembly of independent interlayer dimers. 
Treating the intralayer coupling 



H h 



■h 



J2 {sfs A 



3 1 e>i Oj , 



<h3>r. 



E {#8* 



3 ' °i S 3 



<l,3>r 



, (3.1) 



as a perturbation, we apply the expansion with respect 
to z™.^ 1 using the method of Gclfand, Singh and 
HuseEloElP and GelfandO for various values of a. In 
order to calculate the staggered susceptibility xn and 
collinear susceptibility XC; we a l so a< id the following 
magnetic field terms with wave number Q = (it, it) and 
(it, 0), respectively. 



N 



(3.2) 



and calculate the ground state energy E(h) up to the 
second order in h. Here r,i is the position of the z-th site 
and N is the number of the lattice sites in a layer. The 
susceptibility is given by, 

d 2 E{h) 



X 



dh 2 



(3.3) 



h->0 



where \ stands for xn or Xc according as Q = (it, n) or 
(tt, 0). Using the method of Gelfand,E2P we also calculate 
the expansion series for the single particle excitation en- 
ergies An and Ac at the wave vector (tt,tt) and (tt, 0), 
respectively. 

These quantities are expanded as a power series in z 
and a z as 



00 00 



= E E <*«**(<**)* = E c ^ a > 1 



p=0 g=0 



fc=0 



and 



k 
Ck(a) = y^ c k ~ q , q a q 

q=0 



(3.4) 



(3.5) 



Here O stands for XN: XCj ^N an d A^. Actually, the 
coefficients Ck(a) are calculated up to the 6-th order in 
/3 _1 = J1/J3 for 7 different values of a and c Pi9 's are 
calculated by inverting the relation (p\5|). 



The ratio series of these series are, however, ill-behaved 
except for the close neighbourhood of a — 0. In order to 
locate the phase boundary as precisely as possible from 
the limited data, we assume that the phase transition 
of the present model belongs to the universality class 
of 3-dimcnsional classical Heisenberg model for which 
X ~ (/3 rfc)~ 7 and A ~ (fi - (3 C ) U with 7 ~ 1.4 and 
v ~ 0.71H even in the presence of frustration. This is 
expected to be valid because the Berry phase term always 
cancel between the two layers even if it exists in the single 
layer model and the remaining long wave length action 
is given by the 3-dimensional 0(3) nonlinear a model,cj' 
Thus we obtain the biased [L, M] Pade approximants for 
each value of a as, 



O x ' x [L,M] 



£-,1=0. L Pi 



1 %' 



(3.6) 



with gg = 1 and L + M < 6 where A stands for 7 and 
—v. From the poles z c of the approximants for x~N an d 
Aj^ (xc an d Ac), we determine the critical values (3 C — 
l/z c of the phase transition between the Neel(collinear) 
phase and the spin-gap phase for each value of a. These 
approximants behave as 



1/x [L,M] 



A? 



B? 



P- 



(3.7) 



For 



/3c' 

in the neighbourhood of poles z — z c . Depending on 
L and M, we find many poles which are rather scat- 
tered. Among them, we only accept the positive poles 
with smallest z c (largest (3 C ) and positive amplitudes. 
The poles with amplitudes A° less than the cut-off value 
£a = 0.01 ~ 0.1 are discarded as spurious. Figures 1(a), 
(b) and (c) show the a-dependence of the poles of the 6- 
th, 5-th and 4-th order approximants with L = M, Mil. 
For small a, the critical value of f3 decreases linearly 
with a. This behavior is common for all poles.! 
Fig.l and consistent with other calculations rl 
general values of a, it is physically reasonable to as- 
sume that the critical value of j3 decreases (increases) 
with the increase of a for Neel(collinear)-spin-gap tran- 
sition. Some poles, however, show the opposite behavior 
as shown in Fig. 1. We assume these poles are physically 
meaningless. If these poles are omitted, the qualitative 
features of the phase diagram is common for all approx- 
imants. Namely, no acceptable real positive poles are 
found in the interval 0.45<a<0.65 indicating that the 
spin-gap phase is stable for the single layer model in this 
interval of a. The Neel (collinear) ordered state appears 
for a < 0.45 (a > 0.65). This is consistent with the exact 



diagonaliz 
mations,B' 




and some approximate esti- 
_!l_f although the precise value 
of the critical a depends on the method used. On the 
other hand, the corresponding amplitude B^ does not 
show any singular behavior as j3 c — > on these poles 
as shown in Fig. 2 for the [3, 3] approximants of An 
and Ac- This suggests that the universality class of the 
transition in the single layer model belongs to the same 
universality class as the bilayer model. This is consis- 
tent with the prediction that the Berry phase term is 
dangerously irrelevant even if it exists in the single layer 
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model@) 

Using the [3,3] approximant, the /^-dependence of the 
energy gap is shown in Fig. 3 for various values of a. It is 
clear that the energy gap A^ at (71", tt) increases and A^ 
at (71", 0) decreases with a. The crossover point a cr (/3), at 
which the position of the smallest gap shifts from (ir, ir) 
to (it, 0), varies with (3 as shown in Fig. 4. For large val- 
ues of (3, a cr is close to 0.5 at which the classical ground 
state changes from the Neel state to the collinear state. 
It shifts to 0.576 as (3 tends to 0. It should be noted 
that the phase boundary between the Neel phase and 
the collinear phase is also-siuftcd to 0.6 in the modified 
spin wave approximationOL^P for small (3. This can be 
interpreted in the following way. The dominant short 
range order is Neel type or the collinear type according 
as a ;< 0.576 or > 0.576 for small j3. In the modified spin 
wave approximation, the corresponding long range or- 
der is established because of the underestimation of the 
quantum fluctuation. 

§4. Summary and Discussion 

The spin-1/2 bilayer J1-J2 model is studied by means 
of the dimer expansion method and the ground state 
phase diagram is obtained by the biased Pade analy- 
sis assuming the universality class of the 3-dimensional 
Heisenberg model. For small interlayer coupling, the 
critical value of (3 for the transition between the Neel 
phase and the spin-gap phase decreases linearly with a. 
Within the available data, the spin-gap phase remains 
stable down to j3 = for 0.45< ; a<0.65 which is consis- 
tent with some of earlier estimations. It is also suggested 
that the phase transitions in the bilayer model and the 
single layer model belong to the same universality class. 

The excitation gaps at (tt, it) and {-k, 0) are calculated 
as a function of a and f3 using the [3,3] Pade approx- 
imant. It is shown that the minimum gap shifts from 
(7r,7r) to (ir, 0) at a cv which is close to 0.5 for large (3 
and grows to 0.576 as (3 — > 0. 

At the first glance, these results appear to contradict 
with, the results of the modified spin wave approxima- 
tion!]' Ef which predicts the absence of the spin-gap phase 
in the single layer model. This method also predicts sub- 
stantially large critical value of f3. These are due to the 
underestimation of the quantum fluctuation in the mod- 
ified spin wave approximation. From this point of view, 
the present results are consistent with the modified spin 
wave resultsEiJ if the long range orders found in the lat- 
ter approximation is reinterpreted as the corresponding 
short range orders. 

Needless to say, the present conclusion is far from con- 
clusive. The order of the expansion is still too low and 
only small number of approximants are available. Higher 
order calculation is required to obtain more reliable re- 
sults. Unfortunately, however, the number of the dimer 
expansion graphs becomes enormous due to the presence 
of next nearest interaction. For example, it amounts 
64303 even for k = 6 (present calculation) and the CPU 
time consumed for the calculation is nearly 10 hours on 
FACOM VPP500 supercomputer for each a. For k = 7 
the number of the graphs increase by more than a fac- 
tor of 10 and the calculation of each graph requires even 



more computational time. 

The numerical simulation is performed using the FA- 
COM VPP500 at the Supercomputer Center, Institute 
for Solid State Physics, University of Tokyo and the 
HITAC S820/80 at the Information Processing Center 
of Saitama University. This work is supported by the 
Grant-in-Aid for Scientific Research from the Ministry 
of Education, Science, Sports and Culture. 
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Fig. 1. The poles (3c = 1/zc for (a) N = 6, (b) N = 5 and 
N = 4. The symbols are defined in the figure. The points in the 
left (right) half of the figures are the poles of Xi\j an d A^r 
(X^ 1 and A- 1 '"). 

Fig. 2. The amplitude B^ for [3, 3] approximants of Aj^ and Aq. 
The symbols are common with Fig. 1(a). 



Fig. 3. The /3-dependence of the energy gaps Aj^ and Aq for 
various values of a based on the [3,3] Pade approximant. The 
symbols arc defined in the figure. 

Fig. 4. The crossover point ocr at which the minimum energy 
gap shifts from (n, it) to (w, 0) based on the [3,3] Pade approxi- 
mant. 
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Fig. 1(b) K. Hida 




Fig. 1(c) K. Hida 
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